Anderson Localization in Euclidean Random Matrices 
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We study spectra and localization properties of Euclidean random matrices. The problem is 
approximately mapped onto that of a matrix defined on a random graph. We introduce a powerful 
method to find the density of states and the localization threshold. We solve numerically an exact 
equation for the probability distribution function of the diagonal element of the the resolvent matrix, 
with a population dynamics algorithm, and we show how this can be used to find the localization 
threshold. An application of the method in the context of the Instantaneous Normal Modes of a 
liquid system is given. 

PACS numbers: PACS 



Introduction. Some forty years ago Anderson 0| 
pointed out that disorder can turn a system expected 
to be a metal from band-theory (i.e. a gapless system) 
into an electric insulator. This gave birth to the difficult 
problem of Anderson localization H, H, |^ . The phys- 
ical picture is roughly the following: In a pure system, 
states are described by Bloch wavefunctions. Disorder 
divides the energy band in localized (wavefunctions ex- 
tending over a limited number of unit cells) and extended 
regions (wavefunctions are plane-wave like, involving an 
extensive number of lattice sites). The energy marking 
this division is called mobility edge or localization thresh- 
old. If the Fermi level lies in the localized region, trans- 
port is strongly hampered. Given the large variety of 
systems that are described to some degree of approxi- 
mation by random matrices , it has become clear that 
Anderson localization is significant beyond the physics 
of disordered metals. However, there are very few sit- 
uations 0, where it is possible to estimate the posi- 
tion of the mobility edge in a controlled manner, and 
often heuristic estimates are used (some authors identify 
it with the loffc-Regel limit or with the limits of the spec- 
trum obtained in some kind of effective-medium approx- 
imation 0), or numerical evaluations based on a variety 
of ideas (see e.g. refs.l^orlu?). In this Letter we present 
a method to locate the mobility edge for Euclidean Ran- 
dom Matrices [m (ERM), which we illustrate in a liquid 
calculation of Instantaneous Normal Modes (INM) . 

The entries of an ERM are deterministic functions of 
(random) particle positions. Conservation laws (sec be- 
low) are encoded in the constraint that the sum of all 
elements in a row is null. ERMs appear in the study 
of disordered d-wave superconductors lll^ . disordered 
magnetic semiconductors |l4j (very similar to a spin- 
glass model 0), INM in Hquids vibrations in 
glasses the gelation transition in polymers ^3 ^^'^ 
vibrations in DNA Also, theoretical studies have 
been carried out l2li . 



ERMs in liquids. Consider a translationally invariant 
system with Hamiltonian H[x\ = 'Ylii<j^i.^i ~ "^j)-! fo'^ 
some pair potential v{x). The system has N particles in 
a box of volume V ^ the particle density being p = N jV . 
We impose a long-distance cut-off on w(x), as usual in 
numerical simulations. In studies of short-time dynam- 
ics |l2j the harmonic approximation is made around the 
oscillation centers x'^ . We thus face the Hessian matrix 



where w^i/(af) is the second derivative matrix of v[x). 
As it should, the sum of all elements in a matrix row 
is zero. Once the probability density function (pdf) of 
the oscillation centers, P[a;'^], is specified, the study of 
the spectral prop erties of the Hessian is a problem on 
ERM theory [llj. In the context of supercooled liquids, 
people have considered three types of P[a:'^|: equilibrium 
configurations (the INM), minima |22l | and saddle- 
points [i^l of the potential energy surface. As in all pre- 
vious analytical studies, we restrict ourselves to INM. 

The resolvent and its equation. We need to consider 
the resolvent matrix 0, Ri^,jv[z) = {z — M)~^ {z is 
complex). To study eigenvalues clustered around A, we 
will set z — X + le for small e. Our goal is to find an 
equation for the pdf of the diagonal term of the resolvent 
matrix 7^[Rii(z)], focusing its imaginary part. Con- 
sider the representation of Rii(z) in terms of the eigen- 
vectors la) of the matrix M: 



Im Rjj (A -|- ie) = Im ^ 



I01«)l^ 



A + ie — Aq, 



(1) 



Previous calculations d, 0, 0| only addresed the mean 
value of 7'[Rii(z)], which gives the density of states 
(DOS) through g[\) = Im Tr R,,(A + iO+) (the 

overline indicates mean value). Localization studies [2| 
require at least the variance of 'P[Ri2;(2:)], as can be seen 
from Eq. ^ Consider a A in the localized region. The 
amplitude of eigenvectors |a) with |Aq — A| < e will be 
large for the main particle of the eigenmodc, and will de- 
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crease exponentially with the distance. Let n(A) be the 
typical number of particles for which the amplitude of 
an eigenmode is sizeable. The probability that in Eq. 
particle j will significantly participate in an eigenvector 
of energy Aa G (A — e, A + e) is of order e(7(A)n(A). In this 
case ImRjj is of the order l/(en(A)), while it is of order 
e with probability 1 — sg{X)n(X). Thus, the mean value 
of ImRjj is finite for small e and proportional to g{X), 
but its variance is (t^(A; e) ~ g{X)/{en{X)) and it diverges 
for vanishing e. A more refined analysis 01 shows that 
the pdf for ImRjj decays as {Im Rjj)~^ , (3 < 1.5, with an 
unknown cutoff function that prevents Im > l/e. On 
the other hand, in the extended region, all the | (j |q;) p arc 
0(1 /TV) and one can replace the sum with an integral. 
It follows that the typical value of Im Rjj is a number of 
order one when e = 0"*". 

To get an equation for V[ILii{z)], we have considered 
the generalization of the Cizeau-Bouchaud recursion re- 
lation |7| to the case of translationally invariant sys- 
tems • This equation relates Rii(z) in a system with 
N particles to the full resolvent of a system where particle 
i is not present anymore, but the positions of all other 

— 1 particles are as in the original system. To use 
that equation we need to make a rather drastic simplifi- 
cation regarding particle correlations: following ref. Il^ 
we regard the matrix elements Mi^tjv[x'^] as independent 
random variables, P[M] = Y[i<j P(My )- Thus we throw 
away most of the information contained in the Boltzmann 
weight -P[x'^], restricting ourselves to the pair correlation 
function g'^^^r). Therefore, p(My) is the pdf of the sin- 
gle matrix element, 

T 



P(My) 



1 



P_ 

N 



drg^-'\r)S{M,,+V'ir)) , (2) 



where 7 = inp dr' r''^g^'^\r') is the average num- 
ber of particles whose distance from particle i is less than 
the cut-off. In other words, we have considered the prob- 
lem on a random graph 0, |2^: each particle is the 
root of a Cayley tree whose fluctuating connectivity is 
distributed with a Poisson law of mean value 7. Under 
these hypothesis the central equation of this Letter is 
obtained El: 



R-M ~ (3) 



z + 5]M,,,-^M,,(R^i+M,,)-iM,, 



Here the inverse-matrix symbols refer to the 3x3 dimen- 
sional space of the directions for particle displacements. 
Note also that all the R~-^ are independently choosen 
from the 'P\R.jj{z)] distribution. 



Solving the equation. We shall consider the case of 
the INM of a monoatomic system with potential (in 
natural units) v{r) = (1/r)^^, with a smooth cutoff 
at r™* = 15V^/13 (see Grigera et al. in for de- 
tails). As the only relevant thermodynamic parameter 
is F = pT^^/^, we shall take p = 1. We have obtained 
the pair correlation function by means of a Monte Carlo 
simulation of an iV = 2048 system, for F between 0.3 
and 1.1, in the liquid phase (cristalization happens at 
Fc ~ 1.14). We have also obtained the INM spectrum 
numerically (diagonalizing 100 Hessian matrices), in or- 
der to assess the quality of our approximations for parti- 
cle correlations. Once the g^'^\r') is known, one needs to 
solve Eq. © . We have done this numerically, by means of 
a population dynamics algorithm (PDA). Starting from 
a population of N elements with Rii(2;; t = Q) = I/z, we 
iterate the scheme 



Huiz-t + l) = F(Rjj(z;t),M,j) 



F(R,„M, 



5]M,,(R-/ + My)-'M, 

3¥=i 



(4) 
(5) 



To form the sums in (0, we divide the space around 
particle i up to the cutoff distance, in spherical shells of 
width Ar = 7-'^/4096. The probability of having a particle 
in the sheU is Airp J^^'^W r'^g^'^'> {r'). Both the identity, 
j, of the particle interacting with particle z, and the di- 
rection of the vector r = r ^ — rj are choosen randomly 
with uniform probability. A single time-step consists on 
the sequential update of the full population. We equili- 
brate for 100 time steps, which is verified to be enough 
by checking the time evolution of the first moments of 
the distribution. We then evaluate ■p[Rii(2)] by means 
of (typically) 300 population-dynamics time steps. The 
population averages of Imi?i^.i^i and (Imi?i^.i^)^ are cal- 
culated, then time-averaged. Fr om them wc compute the 
DOS, g{X) = —ImRifj^^ifj^/n, and ImRf^ . An important 
test is to compare the DOS obtained from Eq.Owith the 
DOS obtained by numerical diagonalization of the scram- 
bled matrix M, built by picking its off-diagonal elements 
randomly and without repetition from the off-diagonal 
elements of the true Hessian, and then imposing the con- 
straints of symmetry and translational invariance. This 
kills three-body and higher correlations, making Eq. |2| 
exact. In Fig. ^ we show the DOS at F = 0.6 and 1.1, 
computed both from Eq. Inland from the numerical diag- 
onalization of the Hessian. At F = 0.6 good agreement is 
found. At F = 1.1 we also show the DOS of the scrambled 
Hessian, which agrees with Eq. |3| One can understand 
the mild disagreement between Eq. and the INM spec- 
trum at F = 1.1 by comparing the pdf of the diagonal 
term of the Hessian matrix as computed from Eq. [21 and 
from the actual liquid configurations (Fig. QJ. This pdf 
coincides with the DOS at the leading order in pcrtur- 
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FIG. 1: Top. Left panel: DOS of a soft-sphere system at 
r = 0.6 from the solution of Eq.^lfor e = 1.0 {N = 20000) 
and from the numerical diagonalization of the Hessian matrix. 
Right panel: distribution of the diagonal terms of the Hessian 
matrix, from Eq. |21 and from liquid simulations at F = 0.6. 
Bottom: As in the top part, for F = 1.1. The spectrum of 
the corresponding scrambled matrix (see text) is also shown. 
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FIG. 2: l{X;e) vs. e near the lower mobility edge at F = 1.1 
for (from top to bottom) A = -140, -130, -120, -110, and 
-100. Inset: Scaling plot l{\;e) vs. e/ec{X) for A > -130. 



bation theory (the diagonal term, Mi,;, is much larger 
than the typical ofF-diagonal term M^), so we under- 
stand why at r = 1.1, Eq. 131 overestimates the weight of 
the imaginary frequencies. 

Localization. In principle, one can estabhsh the lo- 
calization threshold from the behavior of the vari- 
ance of the imaginary part of R. Consider ^A; e) = 
e(ha.B?) / (ImR)^. From the remarks above it follows that 
in the localized region and for small e, Z(A;e) tends to a 
constant, while in the extended region the variance is fi- 
nite and l{X;e) is 0{e). Fig. Elshows l{X;e) around the 
lower mobility edge. For the higher values of A, the ex- 
tended, 0{e) regime is reached, but as the threshold is 
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FIG. 3: Evolution of the median (top) and geometrical av- 
erage (bottom) of the population of imaginary parts of the 
resolvent (see text) for A = -140, -150, -160, and -200. 
Data for N = lO'^ . 



approached, the regime sets in for lower and lower values 
of e, making the threshold difficult to estimate. To over- 
come this limitation, we have tried the phenomenological 
scaling l{X,e) = /(e/ec(A)), with f{x) ~ x for a; ~ and 
ec(A) -> for A ^ Ath (Pig.^. The edX) (that we ob- 
tained only for A > —130) can be fitted to a power law 
c(A - A)T, with 7 = 21 and Ath = -171.6. This huge 
exponent is perhaps an indication of an essential singu- 
larity at the mobility edge (cf. the analitical results of 

However, the fit is not very reliable, since the it does 
not reach close enough to the critical point (moving 
nearer requires lowering e, but this is computationally 
expensive since the statistics must grow as l/\/s when 
the variance is growing as 1/e). 

In search of a more reliable way to determine the 
threshold, we have tried a different method, based on 
the observation that in the localized region and in the 
£ ^ limit, ImRii = except for a vanishing fraction 
of particles, suggesting that a population of real resol- 
vents is dynamically stable under Eq. [S] The idea Q is 
to start with an equilibrated population of real resolvents 
(setting z = A), add a small imaginary part and evolve 
the population with eq. El to see whether the imaginary 
parts grow (extended phase) or tend to disappear (local- 
ized phase). We have analyzed in this way the negative 
mobility edge at F = 1.1: We obtained an equilibrated 
population of real resolvents, then added an infinitesimal 
imaginary part, which evolved with the linearized (with 
respect to the imaginay part) form of Eq.[Sl We evolved 
the real part with Eq. |5l as if the population were real. 
Since the pdf of ImRj*' is extremely broad, one must 
be careful in the quantity one chooses to examine. In 
ref.i it was proposed to consider the quotient St of the 
geometric mean of the imaginary parts at times t and 
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FIG. 4: Values of the slope from a linear fit of {log It) > for 
t > 10. The mobility edge is estimated at A = —152(5). Inset: 
A'^-dependence of the slope. 



t — 1. St was assumed to be i-independent, and it was 
expected to be equal to 1 at the mobility edge. However, 
the threshold found with this method was far from the- 
oretical estimates. Indeed, as can be seen in Fig. 13 this 
quantity has extremely large fluctuations, that make it 
practically useless to determine the threshold. Instead, 
we have found that the median (let us call it It) of the 
pdf of ImRj-*^ is much more reliable. However, obtaining 
the thermodynamic limit for It requires N that grows 
fastly with t. For t < 15 we have found no differences be- 
tween N = 10'^ and N = 8x 10^. One would expect that 
log It ~ at, with a less than (greater than) zero in the 
localized (extended) phase, and this is indeed what we 
find at large t (Fig.O. For N = 10^ we have performed 
10 trajectories, then estimated the errors with a jacknife 
method. A linear fit of log Jt, for 1 < t < 15 gives values 
a that locate the mobility edge at A = —152(5) (Fig.Q. 

In summary, we have addresed the problem of Ander- 
son localization in Euclidean Random Matrices, specializ- 
ing to a liquid Instantaneous Normal Modes calculation. 
Under the random-graph approximation for particle cor- 
relations, which keeps only the information of the radial 
distribution function, a closed equation was found for the 
pdf of the diagonal element of the resolvent matrix R^^. 
This equation was solved numerically by means of a pop- 
ulation dynamics algorithm. We have shown that the 
random-graph approximation works rather well at low 
density, but significantly overestimates the weight of un- 
stable modes when the system is close to crystallization. 
We have shown how to locate the mobility edge from the 
mean of the square of the imaginary part of Ri,;(z). We 
have studied numerically the stability of a population of 
real resolvents, improving over the numerical method of 
Abou-Chacra et al. Q , and obtaining an estimate in rea- 



sonable agreement with the value obtained with complex 
resolvents. 
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